Monte Carlo Simulation of Small- World Models with Weak Long-Range Interactions 



O 
O 

(N 

X> 
(N 



X. Zhang and M. A. Novotny 
Dept. of Physics and Astronomy; Center for Computational Sciences; P.O. 
Box 5167; Mississippi State University; Mississippi State, MS 39762-516^ 
(Dated: February 6, 2008) 

Large-scale Monte Carlo simulations, together with scaling, are used to obtain the critical behavior 
of the Hastings long-range model and two corresponding models based on small-world networks. 
These models have combined short- and long-range interactions. For all three models, the critical 
scaling behavior predicted by Hastings [Phys. Rev. Lett. 91, 098701] is verified. 
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Traditionally, studies of materials have concentrated 
on either models based on regular graphs or random 
graphs (such as the models on percolating lattices pj, 
the z model 0, 0, or the model in Ref. Q). Re- 
cently, much more attention has been paid to networks 
that are midway between random graphs and regular 
graphs, for example, scale-free networks The critical 
behavior, as well as the transport properties, are differ- 
ent from that of the same model on a regular graph. 
The small world (SW) graph dj is a type of graph 
that interpolates between a regular graph and a random 
graph. These graphs have been widely used in m any ar- 
eas, for example, in parallel simulations 0, H, HJ HH LUl 
and in studies of the Internet @- Since a particu- 
lar material is related to a given graph with interac- 
tions based on the two-body interaction approximation 
|12|. many models of materials have been studied on 
small world graphs. These models include the Ising 
model, Heisenberg model, and random- walker models 
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Such studies show that these models exhibit mean-field 
behavior. 

It is not easy to analyze small-world models theoreti- 
cally because of the randomness of the long-range bonds 
0. Hastings introduced a long-range model to explain 
the critical phenomena in small-world systems analyti- 
cally [2(j. His model, in which long-range interactions 
are distributed to all the spins in the system, has no 
randomness and is easier to analyze. He claimed that his 
model has the same universal behavior as the small world 
model, with some constrains, and presented a criterion 
for the crossover to the underlying mean-field behavior. 
Thus far, there are no experimental or computer simu- 
lation verifications of Hastings long-range model, as well 
as its comparision with other small- world models. In this 
article, Monte Carlo simulations for the 2-dimensional 
Hastings model and two corresponding SW models with 
Ising spins are presented and the critical behavior is an- 
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II. MODEL AND METHODS 

All three models studied start with a 2-dimensional L x 
L square lattice with periodic boundary conditions and 
an Ising spin with Sj±l located on each site, and a nearest 
neighbor interaction constant J. For the first model, the 
small-world Ising model, randomly chosen pairs of Ising 
spins are connected with a bond, a fixed number (uu) of 
these small-world bonds are added. Once a SW bond is 
added to a pair of spins, no other SW bond is allowed 
to be added to these spins. These SW bonds do not 
change during the simulation, i.e. they are quenched. So 
the total number of Ising spins in the model is N = L 2 , 
and the number of SW bonds is less than =^-. Thus, the 
ratio of the number of small-world bonds to the regular 
nearest neighbor bonds is p = jf^. On average, a spin 
in the system has z — 4 + j% = 4 + 4p nearest neighbor- 
bonds. We suppose the interaction via SW bonds is the 
same strength as that via regular lattice bonds. The Ising 
Hamiltonian of such a system is given by 



H 



-J SiSj J SiSj 



(1) 



SW 



The first sum is over the four nearest neighbor spins on 
the square lattice, and the second sum is over the w SW 
bonds. We use J = 1 in this article. 

Hastings constructed his long-range model by giving 
each spin on the lattice a weak coupling, of order % = -rs, 
to every other spin in the lattice, instead of adding long- 
range bonds with probability p [2fJ . Hence the Hastings' 
Hamiltonian can be written as 



H 



-J ^ SiSj 




N 
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Here M is the magnetization density, 

M=-y Si . 



(2) 



(3) 
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We utilized standard Monte Carlo (MC) simulations [28| . 
A Glauber flip probability with the site for an attempted 
update chosen at random is used. In our implementation, 
two random numbers, J"i,T"2, are generated. The first 
number is used to randomly choose a spin for the spin 
flip attempt, the second is used to determine whether the 
chosen spin should be flipped. If 

" exp-Wfc B T + exp -iWfc B T ' ' 

the chosen spin will be flipped. Here T is the temperature 
and fee is Boltzmann's constant (in our units fcB=l). The 
current energy is E id and E nevr is the energy if the chosen 
spin is flipped. 

We also study another small-world model, one with 
annealed SW bonds. The model is on the square lattice, 
with a long-range interaction the same as the regular 
lattice interaction J. At each Monte Carlo spin flip trial, 
randomly choose a spin and for this one update attempt, 
add with probability p a small- world connection between 
this spin and another randomly chosen spin. Calculate 
the energies -E id and E ncw using the four square-lattice 
nearest neighbor spins and the small-world connection 
(if it is added) , with use of Eq. J3J , flip the chosen spin. 
This Monte Carlo process is similar to the spin-exchange 
process used by Racz et al [2|| |3(j • Since the long-range 
random connection is annealed in this model, we call this 
model the annealed small-world model. 

The Monte Carlo simulation is performed on a par- 
allel computer with each processing element running a 
particular temperature for systems with size less than 
N = 256 2 , while up to four processing elements run- 
ning a particular temperature for the largest system size 
N = 384 2 . The SPRNG |lj random number generator is 
used. A number of quantities are measured, but only the 
order parameter (|M|) and the Binder 4th order cumu- 



lant (U4) are reported here, \M\ = J2f= 

M-) 



E 



JY 



and 1/4 = 1 



The summation index j runs over 



the K different configurations generated in the Monte 
Carlo simulation. Up to 128 processing elements were 
used in our simulations. It took about 40 hours per data 
run for system size N = 256 2 and 70 hours for system 
size N = 384 2 . For each temperature, averages are taken 
using K=10 6 Monte Carlo Steps per Spin (MCSS). 

Hastings developed a general scaling method by using 
renormalization group theory combined with mean-field 
techniques for small world systems. He predicted that 
for small-world systems with small p, the critical point 
shifts away from the Ising critical point, from Eq. (9) in 
M, 



T c -T c = A t p 



1/7 



(5) 



Here 7 is the critical exponent for the susceptibility for 
the local regular lattice system, T c is the critical tem- 
perature of SW system, and T c is the critical tempera- 
ture of the local system. Here the local system is a 2c? 



square Ising lattice, so T c = 2.2691 ••• . The prefactor 
A t = (A+) 1 / 7 , where is the prefactor in the scaling 
of the susceptibility \- It is also predicted that below 
the critical point, T < T c , the magnetization in the limit 
N -> 00 is (from Eq. (10) in HJ) 



\M\ = A M T t 



T 

1 - 



T r — T, 



(G) 



Here Am is a prefactor and (3 is the critical exponent for 
the order parameter for the local regular lattice system. 
In the models studied in this article, the local system is 
a 2d Ising model, so (3 — i and 7 = X. Hence \M\ = 



AyT c ~T for T < T c , and the mean field amplitude 

A oc p -1 diverges as the long-range interaction strength 
p (or the ratio of SW bonds to regular lattice bonds p) 
approaches zero. The mean field region for this scaling 
is small pfj| , with mean field behavior only over a region 
given by 



\T r - T\ 



oc pi oc 



\T, - T r 



(7) 



Note that this type of scaling is related to the coherent 
anomaly method pioneered by Suzuki |32| . Also, this 
scaling can only be expected to be seen in systems with 
a small amount of long-range interactions (p small) . For 
SW systems with a large amount of long-range bonds 
(p = 0.25), even if the long-range interaction is weak 
(0.01J, 0.05 J, 0.1 J, 0.5J), this predicted scaling was not 
observed for the system sizes we studied. 

For small p, as the system size increases, the crossing 
point of Ui approaches the point where U4 is the pre- 
dicted value for infinite size mean field models, u!±^ = 
0.2705 • • • 0, We use the temperature where U4 = 
as the estimate for the critical temperature T c in 



(MF) 



our scaling. 



III. RESULTS AND SCALING 

Simulations were performed for p = 
0.00125, 0.005, 0.01, 0.0625 and 0.125 for the Hast- 
ings model and the annealed small-world model. For 
the 2d SW model, simulations were performed for 
z = 4.005, 4.02, 4.04, 4.254 and 4.5, which correspond 
to the above p values. The results show that the critical 
temperature for the 2d SW model is larger than that 
for the Hastings model for the same p. That is what 
is expected since the spreading of pair-wise long-range 
interactions to every spin in the system weakens the 
long-range effect. The critical temperature for the 
annealed SW model is a little bit smaller than that 
of the Hastings model. Fig. shows the critical 
temperature relation from Eq. JSJ) for all these three 
models, as well as the theoretical prediction. For each 
model, the data falls on a straight line of the Hastings 
predicted slope of . The small deviation from the slope 




FIG. 1: (color online) The critical temperature relation from 
Eq. The critical temperature of the Hastings model is 

much smaller than that of the corresponding 2d SW model 
(upper curves), but very close to that of the annealed SW 
model. The data show that Hastings' prediction about the 
critical temperature works well. 




FIG. 2: (color online) The order parameter \M\ vs. t = 1 — ~- 
for the Hastings long-range model. 

of — for the smaller values of p is due to fluctuations 
and statistical uncertainty. For the Hastings long-range 
model, the error estimate from 16 independent samples 
for p = 0.00125 is a = 0.0030 for N = 128 2 and 
a = 0.0045 for N = 256 2 . The error estimates for the 
other systems is of the same order. 

Fig. (|2J) is the simulation result for the order parameter 
vs. the reduced temperature t, where t — 1 — X, for 

the Hastings model. Only the branches for T < T c are 
shown. It can be seen from the figure, the data curve 
for N — 384 2 and p = 0.125 approaches a slope \ line 
in the region near t < 0.04, which is consistant with the 
estimate from Eq. JJJ (for p = 0.125, and the largest 
system size we studied, N = 384 2 , T c = 2.734). The 
deviation from the slope of h at small t is due to finite size 
effects. This is seen by comparing systems for p — 0.125, 
the smaller lattices (N = 128 2 and 256 2 ) deviate from 



FIG. 3: (color online) The order parameter scaling from 
Eq. (JSJ for the Hastings long-range model. 



the N = 384 2 data at small t. For a fixed value of N, 
the region with the expected slope of \ is also not seen 
for values of p that are too small (as for TV = 384 2 and 
p = 0.01). Consequently, the prediction of Hastings [2(j 
with a slope of \ is only seen for large enough values of 
N and values of p which are small (but not too small for 
a fixed TV). 

Fig. © to Fig. (0 arc scaling results for the order 
parameter \M\ using Eq. © for the Hastings model, the 
annealed SW model, and the SW model. From Eq. J7J, 
the mean field critical region is proportional to T c — T c . 
The smaller the p value, the smaller is T c — T c and the 
narrower is the mean field critical region [2jj. In the 
mean field region, the scaled data points collapse to a line 
of slope |. For the Hastings model, shown in Fig. the 
scaled data points fall on a region close to a line of slope \ 
near t = 1 — £- = 0.01 for all p for larger N, and expands 
beyond t — 0.02 on the upper side and 0.008 on the lower 
side for values of p that are larger (p = 0.0625 and 0.125). 
The annealed SW model and the SW model present the 
same behavior, Fig. Q and Fig. (JSJ). It is expected that, 
for a fixed size system, the mean field region will be wider 
as p increases (but is still small enough so that the long- 
range interaction is weak). The scaling result for the 2d 
SW model illustrates this. In Fig. JSJ, for system size 
N = 384 2 , the mean field region of the largest value 
of p = 0.125 expands much further to smaller t than 
that of p = 0.0625 and p = 0.01. For a fixed p, the 
mean field region is wider for larger size systems. In 
Fig. ©, Fig. Q and Fig. ©, for p = 0.125, as the 
system size increases from N = 128 2 , 256 2 to N = 384 2 , 
the mean field region extends to smaller t, or, to the 
area where the temperature T is very close to the critical 
temperature T c . But for much smaller values of p, for 
example, p — 0.01, the mean field region is hardly seen 
even for N — 384 2 . Simulations for larger size systems 
are needed to penetrate the mean field region as Hastings 
predicted [20] for much smaller p. Unfortunately, we can 
not run larger size systems due to computer limits. 
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FIG. 4: (color online) The order parameter scaling from 
Eq. © for the annealed SW model. 
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model and the SW model. It is hard to say for which 
model Eq. JBJ works better. Although for p = 0.125 and 
N = 384 2 , these models present a reasonable scaling re- 
sult, the scaling result for a smaller p, here p = 0.01, 
is not good and the mean field region is not obvious for 
our system sizes and statistics. It can not be determined 
whether this deviation from the predicted scaling is due 
to fluctuations or to finite size effects or some combina- 
tion. 



IV. DISCUSSION AND CONCLUSIONS 

Starting from the square lattice Ising model, we con- 
structed a small-world Ising model and an annealed 
small-world model that have the same critical behavior 
as the Hastings long-range model. These two Ising mod- 
els have been simulated, as well as the Hastings model. 
Hastings' predictions for the critical temperature and the 
order parameter scaling [2(j are verified. Results show 
the predicted scaling relation in |2(| for the critical point 
T c , Eq. J5J, works well for all three models, Fig . Q. The 
mean field behavior predicted by Hastings [2(j, Eq. J^J, 
is seen in the mean field region, Eq. JJJ , for system sizes 
N = 256 2 , N = 384 2 and not too small p (p = 0.0625 and 
0.125), Fig. ©, Fig. Q and Fig. ©. To see whether 
Hastings' predictions work for much smaller p (for exam- 
ple, p ~ 0.01), simulations for much larger size systems 
(N >> 384 2 ) and much higher statistics would be re- 
quired. 



FIG. 5: (color online) The order parameter scaling from 
Eq. © for the 2d SW model. 



From the above analysis, the mean field region pre- 
dicted by Hastings 20| for the order parameter is seen 
for large enough values of N and values of p which are 
small (but not too small for a fixed N) in all three Ising 
models we studied, the Hastings model, the annealed SW 
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